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Abatiact 

We consider the stability to thermal perturbation of static 
models of coronal loops. The effects of coolr radiatively stable 
material at the loop base are included. The linear stability 
turns out to be sensitive only to the boundary conditions assumed 
on tne velocity at the loop base. The question of the 
appropriate boundary conditions is discussed^ and we conclude 
tnat tne f ree-surtac<“ condition (the pressure perturbation 
vanishes), rather than the rigid-wall (the velocity vanishes) , is 
relevant to the solar case. We find that the static models are 
thermally UQfitAblfif with a growth time of the order of the 
coLunai cooling time. The physical implications of these results 
for the solar corona and transition-region are discussed. 


1 



I ' Introduction 

Tho linear theory for the thermal stability of static 
coronal-loop models has been investigated by a number of authorsr 
but Wien differing conclusions. In his original work, Antiochos 
(1979) concluded that the static models were thermally unstable. 
Similar results were obtained by Hood and Priest (1980). 
However, Chiuderi, Einaudi, and Torricelli-Ciamponi (1981), Craig 
ana NeClymont (1981), and McClymont and Craig (1981a, b,c) have 
found that the models are either stable or that the growth rates 
for instability are too small to be physically significant. The 
reason for this difference in the results of the two sets of 
autnois is in their treatment of the base of the loop models. 
Antiochos (1979) and Hood and Priest (1980) have not included 
cool material, T - 10^ K, in their model, for which the form of 
radiative loss curve (e.g., Raymond, Cox, and Smith 1976) favors 
linear stability (Field 1965). In addition, Antiochos has consi- 
dered only perturbations with a vanishing first-order heat flux 
at the base. Chiuderi, Einaudi, and Torricelli-Ciamponi; and 
Craig and McClymont argue that the growth rates for instability 
are very sensitive to these assumptions so that the models can be 
effectively stabilised by either including some chromospheric 
material or by changing the boundary conditions on the tempera- 
ture perturbation; in particular, the models can be metastable 
(zero growth rate) if the temperature perturbation, instead of 
the heat-flux perturbation, is assumed to vanish at the base. 
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It is clear, therefore, that the proper treatment of the 
base is critical for determining to what extent the models are 
eiciier stable or unstable. In the next section we discuss the 
question of the boundary conditions in detail. 

II. EfiiiuibALifin fflUAlifinfi And finundAX^ CfindlXinnA 

From the form of the perturbation equations, e.g., Antiochos 
( 1979 ), Fiel^' J.965), it is evident that there are four 
independent spatial derivatives in the problem and, hence, only 
four boundary conditions can be specified. This is to be 
expected since the full nonlinear equations also require four 
spatiao. boundary conditions for a unique solution (e.g., 
Richtmyer and Norton 1967). Hence, at each end of the loop two 
conditions must be specified. One condition generally describes 
the thermal propjrties of the base; for example, the loop base 
may act as a thermal bath so that the temperature perturbation 
van .shes there, Tib - 0' or it may act as a thermal insulator so 
tnat the perturbation heat flux vanip;hes there. Fib " 
other condition describes the inertial properties of the base; 
for example, it may act as a rigid wall so that the velocity 
perturbation vanishes, vib ■ 0, or a free surface, so that the 
pressure perturbation, Pib * 0. Depending on the physical situa- 
tion, either of the inertial conditions may be assumed, but 
clearly not both. Assuming both conditions is physically incon- 
sistent; it implies thit the loop base is both rigid and free. 
Also, it is mathematically incorrect since, along with the 
temperature condition, it is equivalent to imposing six indepen- 
dent boundary conditions on a problem that admits only four. 
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However^ previous authors (Antiochos 1979; Chiuderlr 
Einaudi, and Torricelll-Clamponi 1980; and Craig and McC.lymont 
1981; have, In fact, imposed both Inertial boundary conditions on 
the perturbation* in their case, they were able to find 
solutions to the overconstrained problem because they considered 
only the restricted class of models in which gravity is 
neglected; and the loop area and the coronal heating do not have 
any spatial variation. Under these simplifications, the static 
equations are autonomous, and the equilibrium models are 
symmetric about the loop apex. Hence, the solutions to the 
f ir sn-ocder equations, i.e., the normal modes of the loop, are 
either purely symmetric or antisymmetric about the apex. But, 
for the simplified system with no gravity, the only possible 
antisymmetric solution to the force equation is that Pj^ vanishes 
identically. Hence, all the antisymmetric solutions, irrespec> 
tive of their boundary conditions, will also trivially satisfy 
the free-surrace conditions; in particular, all antisymmetric 
modes for the rigid-v/all conditions will satisfy the overcon- 
strained problem. 

Previous authors have used this result to exclude the 
symmetric modes as valid solutions. The argument used is that 
the antisymmetric modes are more physical because all the 
perturbed quantities, and P^ vanish at the boundary. 
But it is evident from the discussion above that this argument is 
spurious since, in general, no mode can satisfy all these condi** 
tions. For a realistic equilibrium model that has no special 
symmetry properties, the modes will be neither symmetric nor 
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antisymmetclc. Insteadr there will be two distinct sets of nodes 
corresponding to the two distinct physical conditions of either 
rigid walls or free surfaces/ and no node will satisfy both. 
Hence/ within each set all the modes are equally acceptable/ or 
unacceptable. 

The question remains as to which set of boundary conditions 
best represents the physical situation in a solar loop. Ideally# 
one would like to place the loop boundary sufficiently deep down 
in the atmosphere that conditions there do not affect the 
stability properties of the corona and transition region. This 
can be dene for the thermal boundary conditions by placing the 
loop base at a sufficiently low temperature. Since both conduc- 
tion and radiation decrease rapidly for T < 6 x 10 K# the 
tnermai time scale at the base can be made to be very long 
compared to the time scales in the corona or transition region 
simply by including sufficiently cool material at the base of the 
static models. We expect that in this case the growth rates and 
tne eigenfunctions will be insensitive to the thermal boundary 
conditions. 

However# it is n:.t possible to select the base point so that 
the growth rates and eigenfunctions are independent of the iner- 
tiai boundary conditions. The reason for this is that the sound 
travel time in the corona# chromosphere and even photosphere is 
very rapid compared to typical thermal time scales# such as the 
coronal cooling time. Hence# we expect that the stability 
properties of the static models will# in general# be sensitive to 
the inertial boundary conditions assumed for the perturbations; 
we show below that this is indeed the case. 
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We are left, therefore, with having to decide which is the 
proper inertial boundary condition at the loop base. We believe 
that the correct condition is the free surface: » 0, and 
that the rigid-wall condition is not applicable to solar loops. 
Assuming the rigid-wall condition implies that at some level in 
the solar atmosphere a true physical boundary exists across which 
no mass transfer can occur. One may think of a rigid wall as 
representing material of infinite inertia and, hence. Immobile. 
For example, in the loop models, the magnetic field is assumed to 
provide infinite inertia to motions perpendicular to the field 
direction, thereby justifying a one-dimensional loop geometry. 
However, in our model the plasma itself has essentially no 
inertia since that is the assumption implicit in dropping the 
acceleration terms. Therefore, by assuming rigid-wall conditions 
at the base, we are postulating that there is a distinct physical 
difference between the plasma in the loop and the plasma below, 
i.e., the inertial properties of the solar atmosphere must change 
abruptly at some level. To our knowledge, there exists no such 
level. 

It is important to note that such an inertial change is not 
equivalent to abrupt density or pressure changes, which do occur 
in the atmosphere. In particular, the chromosphere may be 
thought of as a narrow interface between the low density and 
pressure corona, and the high density and pressure photosphere. 
Due to the density difference, the velocities in the photosphere 
will tend to be very small compared to those in the corona. For 
example, the velocities in the photosphere beneath a coronal hole 
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are nejligible compared to solar wind speeds. But the Irnportant 
q entity .'.s the mass flux, and this is the same in both . >^'*')n8. 
Also, Qi;G t > the pressure ratio, the position of the top of the 
photosph:!rc will be negligibly affected by pressure changes in 
tne coronr., but this does not mean that the top of the 
photosphere acts as a rigid wall with respect to the corona. On 
ute contrary, n:aterial moves freely from the photosphere to the 
r 'en;.. and vice versa. .he top of the photosphere acted as a 
ri’jid wall, the mass flux would necessarily vanish there; 
.wid .iiors would be no steady-state flows such as the solar wind 
or jiivt. n flow into sunspots. 


III. Aod £i9£n£im£il£iQfi 

A. Free-Surface Modes 

In this letter we calculate the growth rates and eigen- 
functions only for the simplified problem that has been discussed 
by previous authors; more realistic models will be considered 
in a subsequent oaper. In addition, we consider only the so- 
called thermally isolated static model, in which the heat flux is 
assumed to vanish at the loop base (Vesecky, Antiochos and 
Underwood 1979). However, we investigate in detail the effects 
of adding cool material to the coronal loop by including a model 
chromosphere, and we include in the analyses the first symmetric 
mode, which we have omitted previously. For comparison we calcu- 
late the growth rates for both the free-surface and the rigid- 
wall-boundary conditions; although as discussed above, we believe 
tnat the free-surface condition is the physically relevant one. 
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Under the assumptions abovor and neglecting the acceleration 
terms (Field 1965), the fir'Jt-order equations can be written in 
L agrangian form as: 


C* + Q(x) 5 + X C - ?! ( .4 XTq + G(x;)/Pq (1) 

where: 

C Tj , (2) 

Q(x) ■ (n^A + Tq 6 Tq - n^A ) 

- no 6/6 no ( e/no))/(10"® , (3) 

G(x) » (£/no - no 6/6 no (e /no> )/ ( 10“^noTo> , (4) 

X » - lo"^ v/Po . . (5) 


In (1) - (5) P is the total plasma pressure; n is the electron 
density; A is the radiative loss coefficient (e.g., Raymond, Cox 
and Smith 1976); e is the coronal heating function; the 
coetficient of thermal conduction is that of Spitzer (1962); the 
subscript "o* refers to the static model; k is Boltzmann's 
constant; v is the growth rate, i.e» all first-order variables 
are assumed to vary as e^^; and the prime indicates differentia- 
tion with respect to the independent variable x, which 
corresponds to the total number of electrons initially in the 
loop and measured from the top (see Antiochos 1976 or McClymont 
and Craig 1981b). 

Equation (1) roust now be solved for the various boundary 
conditions. The thermal boundary conditions are straightforward 
to implement: if the temperature perturbation is assumec to 
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vanish at the baser then ■ 0} if the heat flux perturbation 
vanishesr then * 0. This is the main advantage to using a 

Lagrangian formulation. In Eulerian coordinates the boundary 
condition of a vanishing perturbati i heat flux has an integral 
form in general; whereas it has the simple differential form 
above in Lagrangian coordinates. Note that there is a physical 
difference between the two formulations. In the Eulerian case we 
are setting the boundary at a particular point in the loop# 
whereas in the Lagrangian case we are setting the boundary at a 
particular plasma particle. (Of course^ both formulations are 
identical in tne case of rigid-wall conditions). Since we will 
show below that the stability properties are insensitive to the 
particular position of the boundaryf as long as it is somewhere 
in the chromosphere, we may as well use the Lagrangian 
formulation for mathematical convenience. 

The free-surface condition is also straigntf orward. Since 
is a constant for the case of no gravity, the free-surface 
condition Implies that the righthand side of (1) vanishes. 
Hence, we are lei.r with a standard Sturm-Liouville problem as in 
Antiochos (1979), which we solve in exactly the same manner as 
before, "he rigid-wall condition, on the other hand, is more 
difficult :o implement. The simplest procedure is to express it 
as an integral constraint, equivalent to Equation (2i>) in 
Antiochos (1979). Letting the position of the base points be ±X, 
(±L in Eulerian coordinates), we find that the perturbed 
pressure, can be determined in terms of the eigenfunctions: 


9 



dx 


(6) 


X 

Pi - (k/L) j C T/3/2 
-X 

Theretore, for the rlgid^-wall condition# we must solve the 
inhomogeneous problem (1)# subject to the constraint (6) (c£. 
McClymont and Craig 1981c). Note that thif is not a standard 
Sturm-Liouville problem; hence# none of the standard theorems 
(e.g.# Morse and Feshbach 1953) on the behavior of the 
eigenfunctions and eigenvalues need apply. Since for this prob- 
lem tnere are two parameters# X and Pj^# that must be determined# 
we do not use a "shooting” technique as in the homogeneous case 
to obtain solutions. Instead# we consider (1) and (6) as 
defining a nonlinear problem for X and and use basically 
Newton- Raph son iteration. Convergence was rapid# i 20 
iterations# for all the cases we ran. 

There are two functions in the problem that must be speci- 
fied: the radiative ^oss coefficient# A (T) # and the energy input 
per particle# £ (n#T)/n. For the former we use a smooth analytic 
approximation to the curve derived by Raymond# Cc - and Smith 
(1976) . The stability properties of the models are not sensitive 
to the exact form of A(t) as long as this form is such that it 
implies radiative instability above ' 10^ K and stability below ^ 
10^ K. 

Since the mecbc^nism for coronal heating is not known# the 
form of the energy input e is essentially arbitrary. In our 
previous work we found that the stability properties are also 
insensitive to the exact form for ^ # at least# for the case 
Where this form is a power law dependence on n and T. Hence# we 
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simply assume that in the corona and transition region, i.e., T ^ 
6 X lO^K, the energy input per particle is constant: e /n ■ 

const. However, for the lower temperatures we use a different 
form for the energy input so that there will be a significant 
mass of material at ^chromospheric” temperature, T < 6 x 10^ K. 
The simplest procedure for obtaining a thick chromosphere in the 
static models is to adjust the energy input rate at low tempera- 
tures so that it becomes almost equal to the radiation loss rate 
(e.g., Craig, Robb and Roll© 1981). Therefore, ve use the 
following form for the energy input: 

e/n ■ Cq tanh(x(T)) + nA (1 - tanh(x(T))) (7) 

where 

X(T) ■ (T/TqJj -1+6)® (8) 

and 

Cq, C]^, Tq^, m and 6 are constants. Equations (7) and (8) imply 
that for 6 << !• e/n -> Cq for T > Tq|j and that c/n -> nA for T “ 
Tq^. The amount of material at Tq^ becomes arbitrarily large for 
arbitrarily small 6 t in particular, one finds a chromospheric 
size scale « 1/6 for the value, m ■ 3, of the exponent in 
Equation (8). Of course, this is a highly unrealistic model for 
tne true solar chromosphere. However, it permits us to 
investigate conveniently the effect of adding cool, radiatively 
stable material to the loop base. 

We have analyzed the stability properties of a wide range of 
static models. In Figure 1 we plot the growth rates of the three 
lowest free-surface modes for a series of static models. These 
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Or with ttithar 


modes are the first symmetric modes, • 

^Irbase ■ 0 ^Irbase " antisymmetric mode, 

^irtop ■ "^Irbase " static models all have almost 

identical coronal properties, i.e., at the loop apex the tempera- 
ture, Tcor ■ density, - 5,5 x 10® cm"^, and 

the half-length, L ■ 1.7 x 10® cm. The models were generated by 

2 ' 

fixing: the energy input rate into the corona, ■ 10 

ergs/sec/particle [Equation (7)]; the chromospheric temperature, 
Tq^j ■ ? X 10^ Kj and the exponent m ■ 3. We vary either the base 
temperature (T|^) , or the mass of the base (Tp - and 5 

varies) . 

The growth rates in Figure 1 are expressed in terms of the 
conductive cooling time of the coronal plasma, i.e., we plot V - 
vT^, where is defined as: 


5/2_Po 
10"® T 


cor 


7/2 


(9) 


For the particular values of the loop parameters above, ■ 1.4 
X 10^ sec; however, is Insensitive to the value of g. Hence, 
Figure 1 should be valid for any static model whose parameters 
lie within the observed range of solar values. Note that for all 
the free-surface modes shown in Figure 1, ^ is negative, 
indicating infltdbiilty. 

There are several interesting features of Figure 1. First, 
it is evident that the growth rates of the modes with thermal 
boundary condition, - 0, are insensitive to the value of 
the base temperature or to the mass of chromospheric material. 
In particular, the first symmetric mode has X - -.65; the first 
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antisymmetric mode can be shown analytically to have ^ « 0 (Craig 
and NcClymont 1981). The modes with a vanishing-perturbation 
heat fluXf on the other hand, are sensitive to the value of the 
base temperature, in agreement with our previous results 
(Antiuuhos 1979); but only if the amount of chromospheric 
material is very small, - .1% of the total loop mass. For models 
witn a significant amount of cool material, the growth rates are 
independent of the thermal boundary conditions or of the amount 
of base material. 

We conclude, theretore, that the static models of coronal 
loops are linearly unstable, at least for the class of models 
considered here. The fastest growing mode in all cases is the 
lowest order one; hence, the instability is a global one 
involving the whole coronal loop. This can be seen in Figure 2, 
V .ere we have plotted the eigenfunction for the case with 
tne largest chromosphere, ^ 20% of the total loop mass. Only the 
mode with Tj^ base " 0 shown; however, the mode with T^ base * 
0 is indistinguishable from it except very near the base. The 
modes exhibit a significant amplitude only in the coronal and 
transition region part of the loop. This amplitude peaks at 
approximately the temperature where the radiative loss rate n^ A , 
is maximum ' 8 x 10^ K. Below this temperature the amplitude 
falls orf exponentially. Although the peak occurs in the transi- 
tion region, the coronal contribution to the "energy" in the 
mode, i.e., /T^^ dx, is approximately equal to the transition 
region contribution. This emphasizes the point that the 
instability is a global one and not restricted solely to the 
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transition region 


B. Rlgid'-Wall Nodes 

In Figure 1 A !■ plotted for the first symmetric rigid-wall 
modes and for each of the two thermal boundary conditions. (The 
antisymmetric modes are identical to those of the free-surface 
case as discussed above.) It is evident that the values of ^ are 
very different from the free-surface case. Firstr they ace all 
CfiflitJLyfif indicating stability, except for the cases with very 
lictle chromosphere, < .1%, and T^^' ■ 0. For these cases the 
value of \ is almost identical to the corresponding free-surface 
mooes. However, for the physically relevant cases with some 
chromospheric mass, the values for X are positive. Note also 
tnat for models with significant amounts of chromospheric mass, > 
10%, the values of X decrease, indicating a tendency toward 
instabiiicy as the chromospheric mass increases. 

The situation can be clarified by inspection of the 
eigenfunction Tj^ (x) for the lowest symmetric mode (Figure 2). 
As in the free-surface case, the mode peaks at T ' 6 x 10^ K, but 
note tnat it has a zero crossing, whereas the lowest free-surface 
mode can have none. The amplitude in the chromosphere is large; 
it dominates the contribution to the "energy” in the mode. As 
the chromospheric mass increases, this amplitude increases and 
becomes more constant. Hence the symmetric rigid-wall modes are 
primarily chromospheric perturbations rather than coronal pertur- 
bations as in the free-surface case. Since our model for the 
chromosphere is highly artificial, the physical significance of 
tne rigid-wall modes is questionable. 


14 



DlflfiUfifilfiA 

The key result of this paper is that the static models of 
coronal loops are thermally unstable# even if the models contain 
arbitrary amounts of chromospheric material. This disagrees with 
the conclusions of several authors (Habbal and Rosner 1979# 
Chiuderi# Einaudi and Tor ricelli-Ciamponi 1981# Craig and 
McClymont 1981 and McClymont and Craig 1981a#b#c)# but agrees 
wicn the conclusions of Hood and Priest (1980) and Hragg and 
Priest (1982). The growth rate of the instability is of the 
order of the cooling time of the coronal plasma# independent of 
the particular parameter of the model. We can identify two 
possible physical manifestations of this instability. One is in 
the observations of velocities in the transition region (e.g.# 
Lices et al. 1976). Since the amplitude of the unstable mode 
peaks in the transition region# one might expect the largest 
velocities to be generated there. Another possible manifesta- 
tion of thermal instability is in observation of cool 
condensations in the corona; pr rticularly# quiescent prominences. 
Our results imply that coronal loops are naturally unstable; and# 
hence# condensations may form spontaneously in the corona. 

Of course# in order to investigate these possibilities# the 
nonlinear evolution of the instability must be calculated. 
Clearly# the instability must not saturate at a low level if it 
is to produce sizable velocities in the transition region or cool 
H-alpha condensations in the corona. As discussed in our 
previous paper (Antiochos 1979)# we believe that the instability 
will result in a nonlinear oscillation of the loop about its 
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scatic state. The period of the oscillation should be related to 
the growth rate of the instability, 10^ sec. for typical loop 
parameters. The amplitude of the oscillation, which must be 
calculated by nonlinear simulation, will determine the physical 
importance of the instability. 

Numerical simulations of coronal loops have been performed 
by several authors. These simulations have not exhibited any 
significant evidence for thermal instability. However, in all 
cases rigid-wall boundary conditions were used, and it is clear 
from our results here that these conditions will tend to 
stabilize the models, at least linearly. In addition, numerical 
effects may be providing unphysical stability to the models. We 
intend to consider, in detail, the nonlinear stability of the 
static loop models in a forthcoming paper. 
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£lau£fi I 

Absolute value of ^ for a secies of static models with 
decreasing base temperature and/or increasing chromospheric mass. 
The results for five modes are shown. Solid lines refer to free- 
surface modes, dashed lines refer to rigid-wall modes. 

Z 

The eigenfunction T^ for the lowest symmetric free-surface 
(solid line) and rigid-wall (dashed line) modes. They are 
plotted as a function of log (r) where r is defined as: r ■ 

X 

J'n^ix) X dx . 
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